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X- 

Lattice QCD results [1] suggest that strongly interacting matter at sufficiently 
high temperature undergoes a phase transition (or a crossover) to a deconfined 
quark-gluon plasma (QGP). Despite the difficulties in identifying clear signa- 
tures of a phase transition in ultrarelativistic heavy ion collisions, recent data 
from the experiments at BNL-RHIC clearly point to the observation of a new 
state of matter [2]. 

Depending on the nature of the QCD phase transition, the process of phase 
conversion in the expanding QGP generated in a high-energy heavy-ion colli- 
sion may proceed in a number of different ways. In general, there is a compe- 
tition between the mechanisms of nucleation and spinodal decomposition [3]. 
Results from CERN-SPS and BNL-RHIC feature what has been called sudden 
hadronization [4] or explosive behavior [5] in the hadronization process of the 
expanding QGP and seem to favor a fast (explosive) spinodal decomposition. 
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Recently, possible signatures of this behavior in high-energy nuclear collisions 
were proposed [6]. 

Most theoretical attempts to understand this behavior focused on the rapid 
changes in the effective potential of QCD near the critical temperature such 
as predicted, for instance, by the Polyakov loop model [7], which would be 
followed by a very fast spinodal decomposition process [5]. In this case, one 
observes a nearly instantaneous decay of the Polyakov loop condensate, pro- 
ducing an exploding source of pions as well as large RMS spatial fluctuations 
in the chiral fields. Previous studies using low-energy chiral effective models to 
investigate how nucleation rates compare to the time scale of expansion of the 
plasma also found very likely that most of the system would reach the spinodal 
region and then undergo an explosive phase conversion [8,9,10,11,12,13]. 

Therefore, both experiment and theory seem to suggest that a QGP gener- 
ated in a high-energy heavy ion collision expands and cools down so fast that 
the processes of chiral symmetry breakdown and hadronization can be de- 
scribed within the framework of an effective potential that is quenched into 
the spinodal region, where long-wavelength fluctuations grow with no barrier 
to overcome. This leads to what we will refer to as the explosive spinodal de- 
composition scenario. However, we will argue that even if the system quickly 
reaches this unstable region there is still no guarantee that it will explode. 
To asses the different possibilities one has to study the time evolution of the 
order parameter of the transition after the system is quenched beyond the 
mean-fleld analysis, in order to incorporate the effects from dissipation, and 
ask how they could modify the explosive decomposition picture. In particular, 
dissipation effects have proved to be important in the context of disoriented 
chiral condensate (DCC) formation in heavy ion collisions [14,15,16,17,18]. 

In this Letter we present an exploratory investigation of the effects of dissi- 
pation in the explosive spinodal decomposition scenario of hadron production 
during the QCD transition after a high-energy heavy ion collision. We use a 
Langevin description inspired by microscopic nonequilibrium field theory re- 
sults to perform real-time lattice simulations for the behavior of the inhomo- 
geneous chiral fields. We show that the effects of dissipation can be dramatic 
even for very conservative assumptions. Analytic results for the short-time 
dynamics are also presented and discussed. 

For the sake of simplicity, we consider an infinite system that is quenched 
to the spinodal and then evolves with a fixed effective potential. We beheve 
that the subsequent evolution of the effective potential should not bring deep 
modifications to our general conclusions once the system has already reached 
the unstable spinodal region. On the other hand, effects from the finite size 
of the plasma will most likely play a non-trivial role [19]. Both effects will be 
taken into account in a future publication [20] . 
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In what follows, wc consider the real-time dynamics of chiral symmetry break- 
down of a QGP created in a high-energy heavy ion collision. We assume that 
the system is characterized by a coarse-grained free energy 



(TX 



where Ueff{(l),T) is an effective potential of the Landau- Ginzburg form whose 
coefficients depend on the temperature, and 4>{x, t) is a real scalar field which 
plays the role of an order parameter that is not conserved, such as the chiral 
condensate. To model the mechanism of chiral symmetry breaking found in 
QCD, we adopt the linear cr-model coupled to quarks, whose standard La- 
grangian can be found, for instance, in Ref. [11]. This approach is widely used 
in the literature and its specificities imply no major limitations to our main 
results. 

Quarks can be treated as fast-moving modes and integrated out using a classi- 
cal approximation for the chiral field, yielding the effective potential f/e//(0) T) 
that is shown in Fig. 1. The pion directions play no major role in the process 
of phase conversion we are considering, so we concentrate on the sigma di- 
rection represented by the field (see Ref. [11] for details). In Fig. 1, we 
show Ueff{4>,T) for a few values of temperature. The critical temperature is 
Tc ~ 123 MeV and the spinodal is reached at Tgp ~ 108 MeV. In order to 
simplify the numerics and to compare with analytic results for the short-time 
behavior of the order parameter evolution, we work with a fit of C/e//(0) by 
a polynomial of sixth degree whose coefficients depend on the temperature. 
Such a polynomial fit is an almost perfect representation of Ugff, and it is the 
fit that is displayed in Fig. 1. 

In our analysis, the evolution of the order parameter t) and its approach 
to equilibrium will be dictated by a Langevin equation of the form 

□ + + = (2) 



where F, which can be seen as a response coefficient that defines a time scale 
for the system and encodes the intensity of dissipation, will be taken to be 
a function of temperature only, F = F(T). The function ^{x,t) represents a 
stochastic (noise) force, assumed Gaussian and white, so that {^{x,t)) = 
and (C(f, t)C{x', t')) = 2 rTS{x - x')S{t - t'). 

Eq. (2) could, in principle, be obtained from a microscopic field-theoretic de- 
scription of the real-time noncquilibrium dynamics of the chiral field at finite 
temperature. This procedure was implemented in the case of a A^"^ scalar field 
theory in Refs. [21,15]. The noise and dissipation terms, which originate from 
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Fig. 1. Fit of the effective potential in the sigma direction for different temperatures 
in the region of interest for the phase conversion, Tgp < T < Tc. 

quantum fluctuations, are engendered by either self-interactions of the chiral 
field or coupling to one or more different fields that play the role of a heat 
bath, provided one incorporates higher order terms in the computation of the 
effective equation of motion for (f){x,t). In fact, it is well-known that one has 
to go up to two-loop corrections in order to pick up imaginary parts in the 
self-energy associated with viscosity and dissipation. Self-interactions of the 
chiral field, as well as its interactions with quarks and anti-quarks, justify the 
inclusion of a dissipation and a noise term such as done in Eq. (2) in the 
framework adopted in this paper. The same is true for an approach that also 
includes the dynamics of the Polyakov loop condensate coupled to the chi- 
ral field and, as we will argue, the inclusion of this effect could dramatically 
modify the results of Refs. [5,7]. 

The main physics behind F is the decay of the cj) field and as such one would 
think that F would be significant only in regions where (f) oscillates rapidly. 
However, this needs not be the case in general. The procedure detailed in 
Refs. [21,15] leads to an effective equation of motion that is richer in struc- 
ture and much more complicated than (2). There is no a priori reason for the 
dissipation function F be so simple and the noise be Gaussian and white. In 
general, one obtains a complicated dissipation kernel that simplifies to a multi- 
plicative dissipation term which depends quadratically on the amplitude of the 
field as ?7i(T)0^(5, t)(f){x, t), where rji is determined by imaginary terms of the 
effective action for and depends weakly (logaritmically) on the coupling(s). 
The fluctuation-dissipation theorem implies, then, that the noise term will 
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also contain a multiplicative contribution of the form t), and be in 

general non-Markovian. The white-noise limit is reobtained only in the limit 
of very high temperature. Nevertheless, assuming F to be a linear function of 
the temperature is a reasonable first approximation, as can be seen from the 
results presented in [15,16]. We will comment further on this point later on 
when discussing our results. 

For the sake of simplicity, we adopt the simple approximate form of Eq. (2) 
for a phenomenological description of the dissipative evolution of the expecta- 
tion value of the sigma field. Although assumedly simple, this analysis allows 
for a clear distinction and comparison of the roles played by dissipation and 
the (explosive) spinodal instability in the spinodal decomposition scenario of 
hadron production during the QCD transition after a high-energy heavy ion 
collision. The simple form of Eq. (2) is also convenient for a comparison of 
our numerical results to (linear) analytic estimates in the region of short-time 
evolution to measure the effect of nonlinearities. 

In our numerical simulations we solve Eq. (2) on a cubic space-like lattice 
with 64^ sites under periodic boundary conditions, with a lattice spacing of 
a — 0.91 fm. We use a semi-implicit finite-difference scheme for the time 
evolution and a Fast Fourier Transform for the spatial dependence [22] . Tem- 
perature is fixed to the spinodal value Tgp ~ 108 MeV. We perform several 
runs starting from different random initial configurations around the infiexion 
point of Ueff which happens at (po ^ 0.162 T and then average the results 
from the different initial configurations. For time steps of At — 0.001/T 
the results become independent of the lattice spacing once it is smaller than 
a ~ 1 fm. 

Before presenting the results of the simulations, it is instructive first to analyze 
the short-time behavior of the solution. One can linearize the equation around 
the infiexion point (po substituting (j) hy (j) = (j)' + (po in Eq. (2) and average over 
the noise. For short times, (f)' is small and the cubic and higher-order terms 
in Ueff{4>') can be neglected, so that the equation for 4>' becomes linear. Since 
the equation is linear, the average over the noise can be done formally. An 
analytic form for the short-time solution of the linear equation can be found 
by using the polynomial fit of Ueff, 



6 



Ueff = ^ On 0' 



(3) 



n=0 



which leads to following equation for the average (0') 




(4) 
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where A = 2a2 + Gas^o + 120,4(60 + 200500 + SOag^o- Note that the constant 
term in [/^'^^ does not contribute to the equation Eq. (4) since the average 
over noise of a constant is zero. One can write the solution of Eq. (4) in terms 
of the Fourier transform (0'(k, t)) of (0'(x, t)) as 

{^'{Kt^ 0)) = Cie^i^'')* + C2e-^^('')* , (5) 



where Ci and C2 are integration constants and Ai(k) and A2(k) are the roots 
of the quadratic equation 

x\k) + rx{k) + (y-\A\)^o, (6) 



where we used the fact that A < in our case. Prom this one sees that for 
short wavelengths, such that ^ r^/4 + \A\, we have (complex conjugate) 
imaginary roots and the solution is oscillatory. For long wavelengths, such that 
k^ < r^/4 + \A\, we have real roots and there is an exponential growth of the 
Fourier components. This exponential growth yields the explosive spinodal 
decomposition. 

Two other limits are also instructive. One is the strong dissipation limit of 
large F, such that the first order time derivative dominates over the second 
order one. In this case, the short-time solution to Eq. (4) is given by 

(0'(k,iRiO)) ^Coe-e^'-l^l)*/^ . (7) 



In this case, one sees that short wavelengths with k^ > \A\ are absorbed by 
the system, while those with k^ < 1^41 explode exponentially. Of course, as 
time increases, (f)' increases and the linear equation is not valid anymore and 
the fully nonlinear equation has to be solved. The other interesting limit is 
the one with F = 0, i.e. no dissipation. We will discuss this limit in connection 
with the solution of the full nonlinear equation in the following. 

We show results of simulations for three different values of the dissipation 
coefficient, namely F = 0, 2T and 4T. It can be argued that the response 
coefficient has the form F(T) 2T/6, where 6 is a number of order one to 
first approximation [23]. The cases considered provide a conservative band 
around the value F(r) ^2T to illustrate the effect of dissipation. 

In Fig. 2 we compare the solutions of the full, nonlinear equation with the 
solution of the linearized equation for F/T = 2. Fig. 2 (a) shows that the roll- 
down for the 0(0'^) potential is slower than for the full Ueff potential. This is 
obviously due to the fact that the falloff of full Ueff is steeper that that of the 
0{4>''^) for small 4>'. At short times, smaller than t = 5 fm~^, one sees that both 
solutions are very close to each other. One interesting aspect of the solutions 
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Fig. 2. (a) The full nonlinear (solid line) and the effective potential 

(dashed line), (b) The corresponding average values of (f) in units of its vacuum 
value (pvac as a function of time for F/T = 2. 

shown in Fig. 2 (b) is that the exponential explosion of the linear solution 
happens much later than the explosion of the full nonlinear equation. This 

seems at first sight very counterintuitive since, from the discussion above, one 
would expect an early-time explosion of the solution of linearized solution. 
This does not happen here because the C(0'^) potential is much shallower 
than the full Ueff for small (f)'. 
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Fig. 3. Average value of the chiral field (j) in units of its vacuum value 4>vac as a 
function of time for T/T = 0, 2, 4. 

In Figure 3 we show the average value of in units of its vacuum vahic, (j)vaci 
as a function of time for the three different values of F mentioned above. In 
the absence of dissipation, the solution (dotted line) is obviously oscillatory, 
with an explosive early-time behavior. The effect of dissipation is of retarding 
the exponential growth, as shown by the solid and dashed lines in the Figure. 

The results clearly show that even for a very conservative value of dissipation, 
F = 2T, the effect can be dramatic. For this value of F, dissipation retards 
the time evolution of towards its vacuum value in ~ 100% compared to the 
case with F = 0. The important point to be noted here is that for expansion 
times of the order of 5 fm~^, which in of the order of the time scales for the 
RHIC collisions, there might be not enough time for the onset of the spinodal 
explosion [5]. 

Of course, effects brought about by the expansion of the plasma [5] and by 
its finite size [19], as well as a more realistic treatment of dissipation from the 
microscopic point of view, will bring corrections to this picture. For instance, 
the authors of Ref. [5] consider a Hubble expansion of the system which intro- 
duces a dissipation-hke term to the evolution equation of the form H(d(p/dT), 
where r is the proper time and H — 1/t is the expansion rate. Therefore, 
for a very rapid expansion, corrections due to dissipation, such as discussed 
here, should play a comparatively less important role. Also, dissipation being 
mainly the result of the decay of the 4> field, common wisdom would suggest 
that its effect should be less important at short times, when the field is slowly 
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starting to roll down the potential. However, as discussed earlier, T depends, 
in general, on - in a model it is proportional to 0^ - and it is not a priori 
clear what will be the effect of such a dependence on the short-time evolution 
of the system. These issues will be addressed in a future publication [20]. 
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